program laurentian3
use my_own_da
type(my_taylor) x,p
type(my_taylor) dx_1,dp_1
real(8) L,DL,k1
integer i,n

my_order=1
dx_1=1.d0.mono.1
dp_1=1.d0.mono.2
x=0.d0 +dx_1
p=0.d0 +dp_1

!   qf: QUADRUPOLE, L:=1,K1:=1;
!   qd: QUADRUPOLE, L:=1,K1:=-1;
!   d: DRIFT, l=1;

!sl: LINE=(qf,d,qd);

n=100

L=1.d0
DL=L/n
k1=1.d0

do i=1,n
x=x+DL/2.d0*p
 p=p-k1*dl*x
x=x+DL/2.d0*p
enddo
p=p-x**2
x=x+L*p

k1=-1.d0

do i=1,n
x=x+DL/2.d0*p
 p=p-k1*dl*x
x=x+DL/2.d0*p
enddo
p=p-x**2

write(6,*) " x "
call print(x,6)
write(6,*) " p "
call print(p,6)

write(6,*) " hit return to end "
pause
end program laurentian3
